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Abstract 
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dispersion relation of the Kelvin modes of the vortex line and also study the coupling between the 
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quadrupole modes, both classically as well as quantum mechanically. The quantum mechanical 
solution reveals the possibility of generating nonequilibrium squeezed vortex states by strongly 
driving the quadrupole modes. 
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I. INTRODUCTION 



Quantized vortices play a crucial role in understanding the stability and decay of su- 
perfluid flow. One way to understand this is that the superfluid flow of the Bose-Einstein 
condensate in a simply connected region is irrotational. Therefore, the superfluid can only 
respond to a rotation by the nucleation of singular topological defects i.e. vortices. Fur- 
thermore, the famous Landau urgument shows that a superfluid moving past an object will 
only experience drag if the superfluid velocity is large enough. Again this drag typically 
manifests itself by the nucleation of vortices. 

In the field of ultra cold atomic gases a single vortex in a Bose-Einstein condensate was 
first observed by Matthews et al. jlj and soon after that also by Madison et al. J2|. Since 
then the field has developed rapidly. The studies of a single or a few vortices 0, US 0, Q 
have more recently also been extended by experiments with Bose-Einstein condensates with 
a very large number of vortices [1, |^ [lE, U, 3|- Apart from some of the most recent 



advances, the review article by Fetter and Svidzinsky [13| provides an excellent overview of 
the extensive theoretical and experimental work on vortices in Bose-Einstein condensates. 

While vortices are very important for understanding superfluid properties of liquids and 
gases, the research into Bose-Einstein condensates in optical lattices has provided an im- 
portant alternative for explorin g th e nature of superfluidity. Bosons in an optical lattice 
realize the Bose-Hubbard model [14'] of which the ground-state is, depending on parameters, 
either a superfluid or a Mott-insulator ^l^, [H, EE ES E3- The quantum phase transi- 
tion between these two phases was indeed recently observed in a beautiful experiment by 
Greiner et al. Furthermore, Bose-Einstein condensates in an optical lattice have been 
used to study experimentally diffraction of matter waves [i^l, number squeezing I21I. col- 
lapses and revivals [22], superfluid dynamics 0, i^l, as well as Bloch oscillations |25l. l26l| . 



There is also a number of important theoretical papers on the dynamics and instabilities of 
a Bose-Einstein condensate in an optical lattices 0, H H H H H 13 • 



In this paper, we study the quantum physics of a vortex line in an one-dimensional optical 
lattice by extending on our earlier work '34'|. This problem is interesting, not only because 
it combines the above discussed two different approaches to the study of superfluidity in 
a Bose-Einstein condensate, but also because the quantum properties of the vortex are 
expected to be more pronounced in a system with a smaller number of particles and in 
lower dimensions. Furthermore, a Bose-Einstein condensate in a one- dimensional optical 
lattice has an intriguingly similar layered structure as the high- Tc superconductors and it 
also appears to be a promising way to reach the quantum Hall regime [2^] , for which we need 
to have about one vortex per particle. Using a variational ansatz for the condensate wave 
function in each site of the optical potential, we derive from first principles the quantum 
theory of the vortex line in an one-dimensional optical lattice and obtain after linearization 
the dispersion relation of the Kelvin modes of the vortex line. Physically these modes 
correspond to a sinusoidal wiggling of the vortex line. It turns out that in an optical lattice, 
as compared to the system without a periodic potential, the theoretical treatment of the 
vortex line is considerably simplified and the dispersion relation can be obtained analytically. 

Kelvin modes or kelvons were observed for the first time in a magnetically trapped Bose- 
Einstein condensate by Bretin et al. [s^ and a theoretical understanding of the observed 
phenomena was provided by Mizushima et al. [s^]. In this experiment, the quadrupole 
mode with azimuthal quantum number m = —2 was observed to excite the Kelvin modes. 
Motivated by this experiment we also use our theory to study the coupling between the 
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kelvons and the quadrupole modes in a one- dimensional optical lattice. We find that the 
coupling is described by a squeezing Hamiltonian [38]. This implies that, when the vortex 
line is treated quantum mechanically, the quadrupole mode will generate nonequilibrium 
squeezed states of the vortex line. In Ref. IHj we have discussed also the spontaneous 
squeezing of the vortex that occurs in equilibrium if also interactions between the kelvons 
are taken into account. In this paper, however, we neglect these interactions and discuss 
both the classical as well as the quantum dynamics of the vortex line when the vortex is 
driven by a strongly excited quadrupole mode. Furthermore, we discuss the dynamics of 
vortex squeezing and how to observe the resulting squeezed vortex states experimentally. 

The paper is organized as follows. In Sec. O we present the theoretical foundations of 
our work. In Sec. Illll we discuss our variational approach and derive the Lagrangian for our 
system. We apply the results of this section in Sec. IIVI where we quantize the vortex line 
and find also the dispersion relation of the kelvons. In addition, we derive in Sec. IIVI the 
quantum theory of the quadrupole modes and solve their dispersion relations. In Sec. Mwe 
study the coupling between the quadrupole modes and the kelvons, and in Sees. IVll and lVllI 
we explore, respectively, the classical and quantum dynamics of the vortex line that results 
from this coupling. In Sec IVIIll we discuss the nonequilibrium squeezing of the vortex line. 
We end with a short discussion of our results in Sec. IIXI 



II. ENERGY FUNCTIONAL 

Our starting point is a cigar-shaped Bose-Einstein condensate trapped by the potential 

V{v) = ^{uy + u:y), (1) 

where ujr and Uz are the radial and axial trapping frequencies, respectively, and m is the 
atomic mass. As we assume a cigar-shaped trap, we further have that Uz -C uJt- The 
Bose-Einstein condensate also experiences an one-dimensional optical lattice 

Vb(r) = ^sin2(^), (2) 

where Vq is the lattice depth and A is the wave length of the laser light. We assume 
that the lattice is deep enough so that it dominates over the magnetic trapping potential 
in the z-direction. When this is true and the number of lattice sites is very large, i.e., 
A -C /z = \/'hlmujzi we can in first instance ignore the magnetic trapping potential in the 
2;-direction. 

The lattice potential splits the condensate into A^^ two-dimensional condensates with a 
pancake shape. We take the lattice to be sufficiently deep such that its depth is larger than 
the chemical potential of the two-dimensional condensate. Although we are thus interested 
in a deep lattice, we consider here only the case that there is still full coherence across the 
condensate array. Specifically this means that the lattice potential should not be so deep 
as to induce a Mott-insulator transition. Typically the required lattice depth to reach the 
Mott-insulator transition in a three-dimensional lattice with a filling factor of one is of the 
order of 10i?r, where is the recoil energy of an atom after absorbing a photon from the 
laser beam. In a one- dimensional lattice the number of atoms in each lattice site is typically 
much larger than in a three-dimensional lattice and the transition into the insulating state 
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requires a much deeper lattice [17[. The Mott-insulator transition in a one-dimensional 
optical lattice has very recently been observed but also in that case the filling fraction 
is of the order of one. 

We use trap units from now on, i.e., the unit of energy is hur, the unit of time is 1/uJr, and 
the unit of length is Ir = ^/h/moJ^. The Gross-Pitaevskii energy functional, which describes 
the gas at low temperatures rotating at frequency Q, is then 



V 1 . 2 2\ K)(r) 



T 



2B 



(3) 



where T"^^ is the two-body T-matrix. In the above units the latter is related to the three- 
dimensional s-wave scattering length a through T^^ = Aira/lr- In addition. 



d 



d 



dy ^ dx 



(4) 



is the angular momentum component in the z-direction. 

For a deep lattice potential it is natural to expand the Bose-Einstein condensate wave 
function in terms of Wannier functions that are well localized in the sites. More precisely, 
we expand the Bose-Einstein condensate wave function as 



(5) 



where n labels the lattice sites and z„ = Xn/2lr is the position of the nth site. For now 
we do not specify the wave functions {x, y) of the two-dimensional condensates, but for 
the Wannier function in the z-direction, $ (z), we use the ground-state wave function of the 
harmonic approximation to the lattice potential near the lattice minimum. This harmonic 
trap has the frequency 

97r , 

(6) 



27r r^— 



m 



and the wave function ^{z) is thus given by 



exp 



2/i 



(7) 



where II = \/h/mujL- 

Substituting the above ansatz into the energy functional and ignoring all but the nearest- 
neighbour interactions, we get the energy functional 
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= J dxdy^Y.^:{x,y) 

<n,m> J 

where (n, m) indicates nearest neighbours, and 



\^nix,y)\'^ -nL, 



^nix,y) 
(8) 
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(9) 
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is the effective two-dimensional coupling strength. Moreover, J is the strength of the Joseph- 
son coupling between neighbouring sites and we have 



1 32 ^ V„(:) 



2 dz'^ hujr 



^z + X/2lr). (10) 



With these assumptions J is a time-independent experimentally tunable parameter. Ap- 
proximating the lattice potential near its maximum by an upside-down parabolic potential, 
we can calculate the Gaussian integral, with the result 



J ' ^^^^ 



71^ \ UJr J \l 



2 / \ \ 2 r 2 
/ \ \ -j^-i 
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The energy functional in Eq. (jSI) is now almost two-dimensional. The third dimension 
is visible only in the last term that describes the coupling between neighbouring layers. 
The energy is characterised by two parameters U2D and J, both of which are experimentally 
tunable. The importance of the on-site interaction term proportional to U2D can be enhanced 
by increasing the number of particles in the sites or by making the lattice deeper. Deepening 
the lattice also decreases the strength of the Josephson coupling J and makes the sites more 
independent. It should be noticed that while J is tunable, it is always positive. Physically 
this means that there is always an energetic penalty for having a phase difference between 
sites. 



III. VARIATIONAL ANSATZ 

We intend to study the coupled system of a displaced vortex line and the quadrupole 
modes. Therefore, a variational ansatz for the two-dimensional Bose-Einstein condensate 
wave function $„(x, y, t) must satisfy two obvious physical requirements. Namely, the ansatz 
for the two-dimensional wave functions in each site must predict the vortex precession dy- 
namics and the quadrupole mode frequencies with sufficient accuracy. The ansatz we used 
in an earlier publication ji^l predicts the quadrupole frequencies very well from the non- 
interacting to the strongly-interacting limit, but unfortunately fails to predict in the strongly- 
interacting or Thomas-Fermi limit the correct vortex-precession dynamics. This has to do 
with the fact that our previous ansatz does not correctly account for the size of the vortex 
core in that limit. This problem is not important for the collective modes, but turns out to 
have serious adverse effects in predicting the vortex dynamics. Therefore, we here use the 
ansatz 



(x, y, t) oc exp {x^ + v'^) - [x^ - V^) - ^xy,n{t)xy + ^ tan 

(12) 

where {xn{t),yn{t)) is the position of the vortex in the nth site and the variational parameters 
^n{t) = e'^{t) + ie'^{t) and exy,n{t) = e^y,n(^) + '^^xy,n{'^)^ that describe the quadrupole and 
scissors modes, are complex. In addition, the wave function in each layer is normalized to 
A^. Since our ansatz does not contain a vortex core, the average kinetic energy would diverge 
in this case. Therefore, we introduce a small distance cut-off ^ as an additional variational 
parameter. Together with the condensate size 1 / \/Bq this parameter will be determined by 
minimizing the equilibrium energy functional. 
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A. Equilibrium solution 



In equilibrium the system has a straight vortex line in the center of the condensate. The 
condensate energy per particle is then given by 

- ' ' ^ + BoU + ^-^ + Ue-^^'T [0, B,e] - (13) 



N 2Bq 2 2Bq 



where 



£/ = ^ (f V/^. (14) 



2tI \lrj V ^ 



Furthermore, F [a, z\ is the incomplete gamma function. Eq. (jl3|) is easy to minimize numer- 
ically for both ^ and Bq, but expressions later on simplify when we note that the short disc- 
tance cut-off is given, with realistic parameter values, to a very good accuracy by ^ = \fBo. 
After this simplification we still have to minimize the resulting transcendental equation for 
the remaining parameter Bq. This is a simple numerical task that must be performed only 
once in the beginning of the calculations. 

The equilibrium solution reveals two physically obvious features. First, the Bose- Einstein 
condensate with the vortex is slightly bigger than the Bose-Einstein condensate without the 
vortex. This is due to the centrifugal potential arising from the vortex velocity pattern. 
Second, the size of the vortex core is smaller for larger condensates. Also, it is good to keep 
in mind that with positive scattering lengths < -Bq < 1 and that in the Thomas-Fermi 
limit of large interactions we have Bq 1. The ansatz in Eq. ()12|) provides an accurate 
description of the system when the interaction strength U is larger than about 10. The 
ansatz fails for very weakly-interacting Bose-Einstein condensates, as then the energy in 
Eq. (1131) ceases to have a minimum. In this regime our previous ansatz can however be 
used l4ll. 



B. Second-order expansion of the Lagrangian 

As we are interested in the kelvon and quadrupole excitations we must expand the sys- 
tem's Lagrangian to second order in e, Cxy, Xn, and yn- To second order there is no coupling 
between the quadrupole modes and the vortex displacements. Consequently, we can treat 
them independently. Also, as we are expanding around the equilibrium solution, we can 
safely ignore the zeroth and first-order contributions. The Lagrangian = — E^, with 
a = QoY K, for the quadrupole modes and the kelvons contain also the contributions that 
result from the general time-derivative term in the action for the Bose-Einstein condensate 
given by 



T [$*, $] = ^ 5^ 1 dxdy - , (15) 



as well as the energy contribution Ea that follow from E [^*, \&]. For the quadrupole modes 
we obtain ^ ^ 

/V" ^ ~9R2 i^n^n + ^xy,n^xy,n) (16) 



2B' 

^ n 
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and 



N 



U ^ 1 
2Bn SBn 



~r I 25 5^ j \ "^^2^''* ^n^xy,n) ' g^2 / ^ j^V^n ^mJ V^n ^mJ 



I / / / \ ( " " V 

\ xy,n ^xy,m) ' \^xy,n ^xy,m) 



For the kelvons we get, in the lowest non- vanishing order in Bq, that 



and 



f = f E(i-r[o.B„1+§) (4 + !/?.) 

^ \ / 

H T 2^ [[xn - + (?/„ -ym) \- 



{n,m) 



(17) 



(18) 



(19) 



In the absence of Josephson couphng the equations of motion for the vortex position describe 
the precession of a displaced vortex around the center of the condensate. With realistic 
interaction strengths the vortex-precession frequency predicted by Eqs. (fTHj) and (IT^ i s 
close to the result using a Thomas- Fermi wave function for the disk-shaped condensate |l3| . 



IV. EIGENMODES OF THE SYSTEM 

With the second-order expansion of the Lagrangian available, we now proceed to study 
the eigenmodes of the system. Our variational ansatz in Eq. (fT^ captures three eigenmodes, 
namely kelvons and two quadrupole modes with quantum numbers m = ±2. In the next 
two sections we solve the kelvon and the quadrupole mode dispersion relations. We deal 
with the nearest-neighbour Josephson coupling by transforming to momentum space in the 
same way as in our previous work |41|]. In momentum space the terms due to the Josephson 
coupling are diagonal and are not more complicated than the other contributions in the 
Lagrangian. 



A. Kelvon dispersion relation 

From now on we choose to interpret vortex displacements as quantum mechanical op- 
erators. This involves the important but natural approximation that the field operator 
iljn{x,y,t), which annihilates a boson in the nth layer in position {x,y) at time t, can be 
replaced with {x,y, Xn(t),yn{t)). Physically this implies that the system can be well de- 
scribed by a classical condensate wave function (x, y) and that all the important quantum 
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fluctations are in the vortex position. Using the results in the previous section and trans- 
forming to momentum space using the convention fk = l/\/N^ 'Y2n^~^''^" fny that in 
momentum space the Lagrangian for the vortex hne becomes 



N 



Bo ^ (^VkX-k - XkV^ 



Bo 



Bo BoT[0,B, 



2 



n + j{k)T [0, 5^] 



{xkX-k + my-k) (20) 



where J{k) = J [1 — cos{kd)] and d = A/2. We proceed by defining the vortex positions 
in terms of bosonic creation and annihilation operators as = (^1.^ + Vk)/2^/l7Bo and 
Vk = "^("^Ifc ~ i'k)/'2y/NBo ji^l. In this manner the Lagrangian takes the desired form 

Lie- = ^ (ivlvk - ujk (k) vlvk^ , (21) 

k 

where the kelvon dispersion is given by 

LUK{k) = Bo/2 + n + T[0, B^] {2J{k) - Bo/2) . (22) 

The solution at zero momentum corresponds to the precession frequency of a straight but 
slightly displaced vortex line. 

It is important to realize that due to the Euler dynamics of the vortex, the coordinates Xn 
and iin are canonically conjugate variables and obey the equal-time Heisenberg commutation 
relation 

[X„,yn] = ^. (23) 

Our theory is valid if the vortex line precesses well within the condensate. For the straight 
vortex line this condition implies that the average displacement a/ (a;^) + (y^) should be 
much smaller than \/\fB^. This leads to a limitation for the expectation value of the 
number operator for the kelvons, namely, 



" k 

Consequently, the total number of condensate atoms sets an upper bound for the kelvon 
numbers. 

Our result for the kelvon dispersion relation has some interesting properties. First, the 
initially negative precession frequency can change sign if the momentum k of the kelvon 
is large enough. Moreover, the dispersion is quadratic for small momenta. Both fea- 
tures are due to the attractive interactions between neighbouring pancake vortices, which 
is harmonic for small separations. This contrasts with the behaviour of the vortex line 
in a three-dimensional bulk superfluid, where the kelvon dispersion relation behaves like 
fc^ln(l//c^) and has a logarithmic dependence for small k 0, 0|. Third, the vortex posi- 
tion is "smeared" by the quantum fluctuations. In particular we have for the kelvon vacuum 
state that (x^) = (y^) = l/A^NBo- Therefore, the quantum properties of the vortex become 
more important in a lattice. This is due to the reduced particle number in every site, as 
opposed to the total number of particles, and the spreading out of the condensate wave 
function as the lattice depth is increased. 
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B. Quadrupole mode dispersion relations 



The quadrupole modes are technically somewhat more complicated than the kelvons, 
but their dispersion relations can be derived in a similar way. The quadrupole modes are 
characterized by the 4-component vector (e^ , , ^'^^y^k ' ^'xy,k) ■ writing down equations 
of motion for all 4 variational parameters, we get a set of coupled first-order differential 
equations. The problem of finding the eigenmodes is then reduced to the problem of finding 
the eigenvectors of a 4 x 4 matrix. In our case this can be done analytically. We find four 
eigenvectors 

/ A{k) \ 



Vi 



—I 
iA{k) 
1 



(25) 



V3 



/ -m \ 

i 

iA{k) 
1 



(26) 



V2 = v]", and V4 = Vg. Of these solutions vi and V2 represent the quadrupole mode with 



azimuthal quantum number m 
azimuthal quantum number m - 



2, and V3 and V4 represent the quadrupole mode with 

-2. The fTinction A{k) is given by 



A{k) 



4^2 + WoJ{k) 



4 + 8SoJ(A;) + Bl{l-W - 2V [0, Bl]) 



(27) 



We can now expand any vector (e^ , e'l , e'^y j. , e'^y^k) into this basis of eigenvectors. Since 
(e^ , e'l , e'^y^k ' ^xy,k) ^ vector of hermitian operators, it is sufficient to deal only with the 
complex expansion amplitudes Oi and 03 of the eigenvectors Vi and V3 respectively. When 
we define the quadrupole creation and annihilation operators with 



and 



03 




- i q2,k) 



£2,-^ + ?-2,fc + i q-2,-k - ^ q~2,kj , 



NA{k) 4 

we ultimately find the desired Lagrangian for the quadrupole modes 

m=±2 k 

where the quadrupole-mode frequencies are then given by 



(28) 



(29) 



(30) 



UBn 



Bor[0,B^] 



+ 2Jik) (^(^) + ^ 



± {Bo - 2fl) 



(31) 
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This expression for the quadrupole frequencies generahzes our earher result in Ref. jSJ] to 
the case of a rotating trap with Q ^ 0. 

If we are to consider the quadrupole modes as a small disturbance, the distortion of 
the condensate due to the quadrupole modes should not be too strong. Quantitatively, we 
require that the energy of the quadrupole modes is much smaller than the ground-state 
energy of the Bose-Einstein condensate. In estimating the condensate energy we can in first 
instance ignore the small contribution from the vortex. In this way we obtain an upper limit 
for the expectation values of the quadrupole number operators as 

^ J2^Ql2,kQ±2,k) « VT+2U. (32) 

^ k 

In Fig. ^ we compare the quadrupole dispersion relations in Eq. ()31|) for momentum 
k = and = 0, with the numerical solution of the Bogoliubov-de Gennes equations for 
the two-dimensional Bose-Einstein condensate. As we explained earlier, our ansatz fails for 
very weakly-interacting Bose-Einstein condensates. Nevertheless, Fig. Q demonstrates that 
Eq. predicts the quadrupole-mode frequencies well within a wide range of experimentally 
relevant interaction strengths. 



V. THE COUPLING BETWEEN KELVONS AND QUADRUPOLE MODES 



In the experiment by Bretin et al. 3^| the kelvons were excited by the quadrupole mode 
with a quantum number m = —2. Therefore, we now proceed to apply our theory to 
investigate also the coupling between the kelvons and the quadrupole modes. To obtain 
coupling between the kelvons and the quadrupole modes, we must expand the Lagrangian 
up to third order. At third order we find important contributions from the kinetic energy, 
the angular momentum term, as well as from the time-derivative term. At third order there 
is also a coupling between the kelvons and the quadrupole modes due to the Josephson 
coupling between the sites, but since we assume a small strength for this coupling we can 
safely ignore this. The dominant contributions are respectively given by 



N 



(33) 



^^^^^^^^ ^ r / ,2 _ 2^ X v] 



(34) 



and 



2^(3) 



E 



yVn-^n ~r Xn yn) ~r ^ \XnXn ynVn) 



(35) 



To obtain the above relatively simple expression for T^^\ we expanded up to first order in 
Bq. Contributions of second-order in Bq depend on the details of the vortex core profile and 
are neglected in the following. 
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FIG. 1: Frequencies of the m = ±2 quadrupole modes in units of cj,. for k = and $1 = and as 
a function of interaction strength. The sohd hne is calculated using Eq. and the open circles 
are calculated by solving the appropriate Bogoliubov- de Gennes equations numerically. 



By using our earlier definitions for the kelvon and quadrupole creation- and annihilation 
operators we can express the third-order contributions to the energy as 



lA{k) 



k,k' 



q2,kV-k-k'Vk' - qi-A+k'^U' 



(36) 



and the third order contribution to the time-derivative term as 

j,(3) 



E 

k,k' 



'A{k) 



q-2,kvl+k'^lk' - q-2,-k'^-k-k'Vk' 



+ q2,kVk'V-k-k' - qi~k4+k'^U' 



(37) 



Using these expressions we can write down and solve the classical equation motions, but 
due to the time-derivatives in Eq. the solution of the operator equations of motion is 
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difficult. However, it turns out that within a wide range of experimentally realistic param- 
eter values we can use the rotating-wave approximation to great effect and thus avoid the 
complication introduced by the time derivatives in Eq. (jHTj) . We, therefore, introduce the 
transformation — ^ ■O^.e"*'^-^*-'^''* and q±2,k (i±2,k^~^'^^^^'^^* and ignore the remaining small 
time-dependent term in the transformed T*^^\ In this way we obtain the interaction picture 
Hamiltonians Hm for the coupling between the kelvons and the quadrupole mode with the 
quantum number m as 



k,k' V 

+ qU-k'^-k-k'Vke''-'^"'''^') (38) 



and 




^2 = E V (-^or [0, i?o^] (k)) (g. 

k,k' V * 



„-iS2(k,k')t 



+ i,-k'vlk'vU) e-'''''''''\ (39) 



where 5-2{k,k') = uj^2{k')—ujK{k + k')—ujK{k) and 52{k,k') = —uj^2{,k')—ujK{k + k')—ujK{k). 
The usual experimental preparation will typically imply that only the quadrupole mode 
with momentum fc' = is excited. Therefore, we assume from now on that only the fc' = 
quadrupole mode is present and drop the momentum index from the operator q2-k'- 

As we show later in Sees. IVII and IVlll the time evolution of the kelvons is strongly 
dominated by the almost resonant modes inside a narrow region in momentum space. Within 
this region the rotating-wave approximation is clearly appropriate. Furthermore, we have 
solved the classical equations of motion numerically without using any approximations. In 
this way we verified that the solution given by the rotating-wave approximation is indeed 
very accurate. 

The coupling between the kelvons and quadrupole modes can only be important if the 
detunings ^±2(^,0) are small. Moreover, the coupling is resonant, when the detuning van- 
ishes. Since in the absence of rotation uoxi'^) is small and negative and ujxik) increases with 
increasing k, whereas uj±2{0) is positive, such a resonance is only possible between m = —2 
quadrupole mode and the kelvons. We demonstrate this graphically in Fig. |21 This reso- 
nance condition singles out the wavelength 27T/kQ for the ensuing wiggles in the vortex line. 
The resonance is possible only if the kelvon frequency increases quickly enough with increas- 
ing momentum k. If, for example, the strength of the Josephson coupling J is too small, the 
resonance is not possible at any value of momentum. In order to have a resonant coupling 
between the kelvons and the quadrupole mode with m = —2 the approximate condition 

^r[o,i?o^]>^ (40) 

must be satisfied. In the absence of rotation the detuning for the coupling between kelvons 
and m = 2 quadrupole modes is always large and therefore the rotating-wave approximation 
eliminates it. 
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For future convenience we define the coupling strength between the kelvons and the 
quadrupole mode with m = — 2 as 



Ec{k, k') = ^^ {BoT [0, B^] -n + UK{k)) . (41) 
In terms of this quantity the interaction picture Hamiltonian becomes 

k,k' 

This interaction picture Hamiltonian is similar to the Hamiltonian encountered in studies of 
parametric processes in quantum optics 3^|. Theories describing the coupling be tween the 



kto,mc aud molecular condensate have al~a similar structure BSlliEliiia. Namely, 

two atoms (kelvons) are combined to form a single molecule (quadrupole mode). 

While the coupling between the kelvons and m = —2 quadrupole modes is more prevalent, 
the situation is slightly more involved in the rotating case. In a rotating trap the frequency of 
the quadrupole modes is shifted by ^2Q. For positive rotation frequencies this shift implies 
that the resonance with quadrupole modes with m = —2 and the kelvons is shifted to a higher 
value of momentum or possibly disappears altogether. On the other hand, the frequency 
of the quadrupole mode with m = 2 is lowered. This signifies that with sufficiently large 
rotation frequencies the quadrupole mode with m = 2 can become resonant with the kelvons. 
In the Thomas- Fermi limit this implies rotation frequencies larger than about 1 / y^. At such 
high rotation frequencies the dynamical instability toward the entry of additional vortices 
is likely to play a role 0, E^- This unaccounted physical process makes the possibility of a 
strong coupling between the kelvons and m = 2 quadrupole modes a technical peculiarity of 
our variational ansatz and renders it, most likely, experimentally irrelevant. 

In this paper we focus on the Bose-Einstein condensate in a one- dimensional lattice, but 
it is clear that the results in this section can also be applied to study an elongated three- 
dimensional Bose-Einstein condensate. The theory for the coupling between kelvons and 
quadrupole modes in such a system is formally similar to the one we presented here. When 
only quadrupole modes with k' = are relevant, there are only two modifications. First, 
the kelvon dispersion relation should be replaced with the correct dispersion relation for a 
vortex line in an infinite cylinder with a radial density profile. Second, the quadrupole-mode 
frequencies should be replaced with the correct quadrupole frequencies for the elongated 
Bose-Einstein condensate. For extremely elongated condensates these frequencies have the 
same functional form as Eq. (j31|) . but rather than using the interaction strength U scaled 
into two dimensions, we must use the correct three dimensional interaction strength. 



VI. CLASSICAL DYNAMICS OF THE VORTEX LINE 

In this section we apply the results in the previous sections to study the classical nonlinear 
dynamics of the coupled kelvons and quadrupole modes with m = —2. The equations of 
motion for the kelvons and the quadrupole mode with m = —2 are given by 

z^ = E,(fc,0)g_2e-^^-('='°V„ (43) 
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FIG. 2: Demonstration of the resonance condition for the couphng between the kelvons and the 
quadrupole mode with m = — 2 when J = 0.1, U = 100, and = 0. The couphng becomes 
important when the kelvon frequency ijJxik) = a;_2(0)/2, i.e., when the hnes cross. Since the 
lowest line does not cross the kelvon dispersion, the couphng between kelvons and the quadrupole 
mode with m = 2 is very weak. The unit of frequency is Ur- 



= -E,{k, 0)qUe^'-^^'''^%, (44) 

and 

= E 0)(4^-.)e^'-('='°)*. (45) 

k 

In the above equations we assumed, an assumption we make throughout this paper, that 
the quadrupole mode is strongly excited and can thus be described with a complex number. 
In this section we also treat the kelvons classically, and thus the operator properties of the 
equations of motion will be, for the time being, ignored. 

We have shown earlier how the coupling between the kelvons and the quadrupole mode 
with m = —2 becomes important at a certain value of the kelvon momentum k. From the 
knowledge of the kelvon and quadrupole mode dispersions we can solve the position of the 
resonance, but we still do not know the width of this resonance. This width is important 
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as it influences how many kelvon modes become occupied. We obtain a simple estimate 
for the width of the resonance by assuming a constant quadrupole mode amplitude ^'-2(0) 
and by seeking exponential solutions for the kelvon amplitudes Vkit) oc exp {Qci,kt). The 
equations of motion for the kelvon amplitudes then result in an eigenvalue problem for the 
characteristic (imaginary) frequency ^lD,k 




This matrix has eigenvalues 

nD,k = ±VEcik,0)q_2iO) [E,{k, 0)g_2(0) ± |5-2(A;,0)|] (47) 

and therefore the Kelvin modes are dynamically unstable only when |(5_2(/c,0)| < 
\Ec{k,0)q-2{0)\. This estimate for the width of the resonance is in satisfactory agreement 
with the numerical solution of the nonlinear time evolution following from Eqs. (j43j) - (j45|) . 
In addition, Eq. (j47|) indicates that the timescale for the onset of the instability is given by 
l/|E,(A;,0)g_2(0)|. 

The argument above was based on the assumption of having a constant quadrupole 
field. This assumption is bound to fail during the time evolution. When this happens, the 
numerical solution of the problem becomes a necessity. In Fig. El we show two examples 
of the classical time evolution of the coupled system. In this figure only the quadrupole 
mode is initially strongly populated. Furthermore, for numerical reasons, the kelvons have 
a small initial amplitude which acts as a seed for the time evolution. When the kelvons are 
treated quantum mechanically such a seed is not required as we show in the next section. 
As is clear from the figure, the quadrupole mode is strongly damped while kelvon modes 
are being populated. At longer times the system exhibits oscillatory behaviour, but such 
oscillations become weaker as the number of lattice sites is increased. 

When the kelvon amplitudes are known, it is a simple matter to calculate also the dy- 
namics of the vortex line. In Fig. |3]we show an example of the vortex line dynamics for the 
system solved in Fig. El (a)- (c). In the figure the short wavelength corresponds to the reso- 
nant mode with the momentum ko and the envelope with the longer wavelength is a result 
of the width of the resonance. The fact that many kelvon modes with different momenta 
are occupied is also apparent from the fairly complex structure of the vortex line. 

VII. QUANTUM DYNAMICS OF THE VORTEX LINE 

In this section we discuss the quantum dynamics of the vortex line. In our numerical 
solution we use again a classical quadrupole mode, like in the previous section, but the 
kelvons are treated quantum mechanically. When the time interval At is short enough so 
that the quadrupole field q-2{t) = \q-2{'t)\e^^^^^^ can be considered as constant, we solve the 
Heisenberg equations of motion for the kelvons from time t tot + At analytically. The result 
is given by 

Vk{t + At) = Vkit) cosh {VL^{t)At) - ivl^it) sinh {Sl^{t)At) e-'W*)+'5-2('='0)t) (43) 

and 

vl^{t + At) = vl^{t) cosh (fie(t) At) + Wk{t) sinh {VL^{t)At) e^C'/'W+^-^C^-o)*), (49) 
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FIG. 3: (a)-(c) The time evolution in units of l/w^ of the total kelvon number per site, the 
quadrupole density per site, and the number distribution of kelvons in the end of the simulation, 
when kelvons are treated classically and Ng = 101. (d)-(f) These plots used the same parameters 
as the plots (a)-(c), except that the number of sites was larger, namely Ns = 751. In both sets of 
simulations we had initially |g-2p/A^s = VT+W/IO and a very small kelvon population to act as 
a seed. Also, we used J = 0.1 and U = 100. In plots (c) and (f) the +-signs indicate the actual 
modes used in the numerical solution whereas the solid line is an interpolation between them. 

where Qcit) = Ec{k,0)\q-2{t)\. This solution enables us to update the correlation functions 
needed later on. In particular, we need the correlation functions {vkV-k) and = {v\vk) 
the first of which is used to update the amplitude of the quadrupole mode. 

In Figs. Elandiniwe show an example of the quantum time evolution of a system that is ini- 
tially prepared in a kelvon vacuum and with a nonvanishing quadrupole-mode amplitude. It 
is instructive to compare Fig. Elwith its classical counterpart in Fig. El (d)-(e). The ambigu- 
ity of the classical initial state influences the growth rate of the classical kelvons. Therefore, 
detailed comparison with the quantum result cannot be easily made. Nevertheless, it is 
clear that in the quantum mechanical simulation the oscillatory behaviour is weaker. Fig. IHl 
shows the kelvon-number distribution as a function of time. After a short transient the 
system settles into a nearly steady state. In principle, revivals of the quadrupole field are 
possible, but the revival time with the number of sites we used in this simulation is longer 
than the length of our simulation. With fewer sites the oscillatory behaviour becomes more 
pronounced, in agreement with the our earlier results of the classical simulations. 
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t=0 



t=20 



t=30 




FIG. 4: Classical vortex line dynamics when J = 0.1, U = 100, and Ng = 101. The initial state 
has \q-2\'^/Ns = ^/T+'2U /lO and a very small kelvon population to act as a seed. The unit of time 
is 1/ur and the unit of length is ^yh/mLOr■ 



(a) (b) 




t t 

FIG. 5: (a) This plot shows the time evolution in units of l/ojr of the total average kelvon number 
per site for quantum mechanical kelvons. (b) This plot shows the quadrTipole density per site as a 
function of time. In both plots we used J = 0.1, U = 100, and Ns = 751. Initially the system was 
prepared in the kelvon vacuum and with a quadrupole density |g_2p/A^s = VT+'2U /lO. 
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FIG. 6: (Color online) This figure shows the kelvon-number distribution as a function of time 
in units of l/w^ when the kelvons are treated fully quantum mechanically. We used parameters 
J = 0.1, U = 100, and Ng = 751. Initially the system was prepared in the kelvon vacuum and with 
a quadrupole density |g__2p/-^s = \/l + 2U /lO. The unit of time is l/tOr- 

A. Transition rate into kelvons using Fermi's golden rule 

Fermi's golden rule 



can be used to estimate the transition rate from the quadrupole mode into kelvons. In this 
expression is the initial state, is the final state and Eq and Ef are their energies, 
respectively. In our case the initial state is a kelvon vacuum with a coherent quadrupole 
mode of amplitude g_2(0), i.e., = |g_2(0), rikg = 0, n^ko = 0). By calculating the 



where ko is the resonant value of the kelvon momentum and = A/2 is the lattice spacing. 

Fig. [7| compares the Fermi's golden rule prediction with the numerical solution of the 
quantum multimode problem by presenting the half-life of the quadrupole field as a function 
of initial quadrupole field strength. As can be seen from this figure, Fermi's golden rule 
underestimates the decay rate for narrow resonances (small values of the initial quadrupole 
field) while overestimating the decay rate for wide resonances. The exact solution deviates 
from the Fermi's golden rule result, since Fermi's golden rule ignores the finite lifetime of 
the final states. In particular, the detunings in Eq. (j45j) contribute to the phase factors 
which determine whether the quadrupole mode is increasing or decreasing. At small times 
the quadrupole mode is always depleted, but at longer times the behaviour depends on the 
width of the resonance. 

VIII. NONEQUILIBRIUM SQUEEZING 

As the interaction between the quadrupole mode with m = —2 and the kelvons closely 
resembles the squeezing Hamiltonians in quantum optics, it is worth while to investigate the 
squeezing of the vortex in detail. This will be done in the following two subsections. We 
start in Sec. IVlll Al by discussing the squeezing when the quadrupole mode is resonantly 
coupled to only one Kelvin mode. Understanding this important special case enables us to 
explore, in Sec. IVlll Bl the experimentally more realistic problem of a multimode squeezing. 

A. Single-mode squeezing 

In the simplest case we have only one resonant kelvon mode at momentum ko. For 
simplicity we assume that initially the quadrupole mode is a real number. We can then 




(50) 



/ 




(51) 
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|9-2(0)|VW. 

FIG. 7: This figure shows the time it takes for the quadrupole field density to halve from its 
initial value |g_2(0)p. The solid line is based on the numerical solution of the quantum multimode 
problem, while the dashed line is based on Fermi's golden rule. We used parameters J = 0.1, 
U = 100, Ng = 751, and the unit of time is l/cOr- 



define two quadrature operators 



1 



and 



2V2 

1 

2V2 



These operators obey an uncertainty relation 



and their commutator is 



bi,b2 



I 

2' 



(52) 
(53) 

(54) 
(55) 
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FIG. 8: (a) The time evolution in units of l/u;^ of the quadrature operator uncertainty y {bf)- We 
used parameters A^^ = 100, N = 100, J = 0.1, U = 100, and |g_2(0)|VA^s = \/l + 2C//10. The 
solid line is for the operator bi and the dashed line for the operator 62. (b) The figure shows the 



time evolution for the product w {byy (62)- 



These operators are analogous to the squeezing operators familiar from quantum optics j38| , 
but are somewhat counter intuitive. We gain more insight when we realize that the operators 
bi are related to the position operators {0) and {0) in a coordinate system rotated by 
an angle 6. In particular, we have 61 = y/WB^Xk^iO) and 62 = ^/NBq ykoiO). The rotation 
angle 9 of the vortex position operators depends on the initial phase of the quadrupole mode, 
because Eqs. and (j^Hj) are in principle modified if ^-2(0) is not real. If the quadrupole 
mode is initially |g_2(0)|e*'^'3 the rotation angle is given by 

^ = (-f + ^q) /2. (56) 

The proportionality of the quadrature operators to the vortex position operators implies 
that the possible squeezing of the vortex state is reflected in the corresponding deformation 
of the vortex-position distribution. 

In Fig. IHl we show an example of the time evolution of the squares of the quadrature 

operators as well as the product \J (^\)\J (bf). As can be seen from this figure, the mini- 
mum uncertainty state (the vacuum) remains a minimum uncertainty state, but can become 
strongly squeezed during the time evolution. In coordinate space this squeezing is reflected 
in the deformation of an initially circular uncertainty ellipse into an ellipse with dramatically 
different main axes. Furthermore, since the initial quadrupole mode amplitude was choosen 
as a real number, the main axes are always along the lines y = ±x. 
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B. Multimode squeezing 



The main difference between the multimode squeezing, as opposed to the resonant single- 
mode squeezing, is that the area of the uncertainty ellipse is no longer constant. This means 
that when the time evolution starts from the kelvon vacuum, the state evolves away from 
the minimum uncertainty state. 

Generalizing our earlier results for the single-mode squeezing we define the multimode 
quadrature operators in terms of single-mode quadrature operators for kelvons with momen- 
tum k, bi^k, as 



|. However, rather than embarking 



For these operators the commutator is again ^i, 62 

on this more formal road, we choose to focus on the shape of the vortex-position distribution. 
This distribution is measurable and provides the necessary signatures of the squeezing. 

When the initial state of the vortex line is a kelvon vacuum, the fluctuations of the vortex 
position in a coordinate system rotated by an angle 9 are given by 



and 



^^^^'^^'^ = ^]V^[l + 2no(t) + 2|(t;o(tH(t))|sin(2^ + 0o(t)) 
+ 5^ (^fc (^) + ^ + I {^k {t)v-k (t) ) I sin {29 + (t) ) 



^^^^'^^'^ = ^^^^[l + 2noW-2|(^;o(iHW)|sin(2e + 0oW) 
+ I] Ufc(^) + ^ - I {vk{t)v^k{t))\ sin {29 + ^t)) 

k^O ^ 



(58) 



(59) 



where we defined {vk{t)v^k{t)) — |('J^jfc(^)f-fc(i))|e*'^''^*^- For this initial state the quantities 
{x{9,t)^) and {y{9,t)^) are independent of the site. Furthermore, their maximum and min- 
imum with respect to the rotation angle define the lengths o'^(t) = Taax[{x{9,ty)] and 
o'y{t) = min [{y{9, t)^)] of the main axes of the uncertainty ellipse as well as their directions 
9max{t) and 9max{t) + 7r/2. In terms of these quantities and in a coordinate system rotated 
by an angle 9max{'t), the uncertainty ellipse is defined by 

+^ = 1- (60) 



The main axes of the uncertainty ellipse are in general time dependent. However, we can 
define the deformation parameter e{t) of the uncertainty ellipse by changing the rotation 
angle 9{t) of the coordinate system as a funtion of time. In this way we define the deformation 
parameter as 

"lit) - al(t) 

^ om^y '''' 
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FIG. 9: This figure shows the deformation parameter of the uncertainty ehipse with two different 
number of lattice sites. The solid line was calculated with Ng = 101 and the dashed line with 
Ns = 751. In addition, we used parameters J = 0.1, U = 100, and |g_2(0)|VA^s = Vl + 2U/W. 
The unit of time is 1 /ur ■ 



In Fig. El we plot the time evolution of the deformation parameter with two different 
number of lattice sites. This figure demonstrates how the uncertainty ellipse first becomes 
strongly deformed and then settles into a more symmetric configuration. In the simulation 
with fewer number of lattice sites the strong deformation of the uncertainty ellipse revives 
at about t = 60. With larger number of lattice sites the behaviour is, at this timescale, more 
irreversible. In Figs. ^1 and ^2 we show snapshots of the uncertainty ellipses corresponding 
to the simulations in Fig. El In these figures we also indicate the main axes of the ellipses. 
The direction of the main axes depend on the initial phase of the quadrupole mode. 



C. Experimental observation 

As noted in the previous subsection the fluctuations {x {0)^) and {y {0)^) are independent 
of the site. Therefore, in each site the fluctuations of the vortex position are described by 
the same uncertainty ellipse and the measurement of the vortex position in each site sample 
the same distribution. However, the number of samples must be large enough to distinguish 
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t=5.0 




FIG. 10: This figure shows snapshots of the time evolution in units of l/w^ of the uncertainty 
ellipse when A''^ = 101. In addition, we used parameters J = 0.1, U = 100, and |g_2(0)p/iVs = 
\/l + 2C/ /lO. Dashed lines indicate the directions of the main axes of the uncertainty ellipse. The 
unit of time is 1 /ujr and the unit of length is ^JhjmijJr- 

the squeezing from statistical fluctuations. The in situ imagining of the three-dimensional 
vortex line is difficult, but not impossible fsi', 1^. 

In many experiments the trap potential is turned off and the system is allowed to ex- 
pand. This makes the experimental imaging of the condensate as well as the vortex easier. 
Therefore, it is important to study how the vortex-position distribution evolves under the 
expansion [sil . Is^ js^l ■ We assume that only the magnetic trap is turned off and that the 
optical potential is left unchanged. When the small inter-site coupling is ignored, each site 
will expand independently in two- dimensions. 

We study the expansion using the same variational ansatz we have used until now except 
that the condensate size parameter Bq is replaced with a time-dependent complex variational 
parameter Bq B'f){t) + iB^lt). The vortex core size changes as the condensate expands 
and we take this into account by using a time-dependent short distance cut-off 

This expression corresponds to the aproximation that the vortex-core size is always equal 
to the coherence length in the center of the condensate. At t = 0, before the condensate 
expands, Eq. (j62|) agrees with our earlier equilibrium value for the short distance cut-off in 
the strongly-interacting Thomas- Fermi limit because Bq(0) ^ l/\/2U. 
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t=0.0 



t=1.0 



t=2.5 





FIG. 11: This figure shows snapshots of the time evolution of the uncertainty ehipse with Ng = 751. 
In addition, we used parameters J = 0.1, U = 100, and \q-2{0)\^/Ns = ^/TT2U/lO. Dashed lines 
indicate the directions of the main axes of the uncertainty ellipse. The unit of time is l/ur and 
the unit of length is ^/Ti/mujr- 

Up to second order in BQ{t) The equations of motion for the variational parameters are 

Xn{t) = C{t)Xn{t) + UJo{t)yn{t), (63) 



(64) 



8C/2 



Kit), 



and 



1 - 



7 - 1 - In 2 - 2C/ - In [/ + 



1 



2U 16U^ 



(65) 



(66) 



To obtain the last equation wc also expanded up to second order in Furthermore, 
7 pa 0.577216 is the Euler-Mascheroni constant, 



C{t) = -B'^it) 



1 

81P 



(67) 
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and 



(68) 



Importantly, from these equations we see that the dynamics of the expansion is independent 
of the vortex position. The real functions C{t) and ujQ{t) depend on -Bo(^) ^"^^ -^o(^)' 
respectively, but since the expansion is independent of the vortex dynamics they can be 
considered as, in principle, given functions. The equations of motion for the vortex positions 
are then solved by 



Xnit) = exp I / C{t)dt 



a;„(0)cos / uJo{t)dt]+yn{Gi)sm\ / iUo{t)dt 



(69) 



and 



yn{t) = exp 



C{t)dt] 
/ 



2/n(0) COS 



uJo{t)dtj 



Xn{0) sin 



uJo(t)dt 1 
J 



(70) 



This result shows that the vortex dynamics during the condensate expansion is a combination 
of precession and scaling. The precession angle and the scaling factor only depend on the 
dynamics of the condensate expansion and are the same irrespective of where the vortex 
was initially located. Therefore, a squeezed vortex distribution expands, to a very good 
accuracy, homologously, i.e., without changing its aspect ratio. The squeezed distribution is 
therefore observable even after expansion of the condensate. 

In Fig. ^1 we demonstrate how the expansion changes the squeezed vortex-position dis- 
tribution by solving the equations of motion explicitly. This figure makes the homologuous 
character of the expansion clearly apparent. 



IX. SUMMARY AND CONCLUSIONS 



We derived the quantum theory of the vortex line in a one- dimensional optical lattice and 
obtained the dispersion relation of the Kelvin modes. Our variational approach enabled us 
to also include quadrupole modes into the system and study their coupling with the kelvons. 
The coupling turned out to be described by the squeezing Hamiltonian and we solved the 
dynamics of the vortex line coupled with the quadrupole mode both classically as well as 
quantum mechanically. Furthermore, we found that the vortex line can indeed be squeezed 
in this system if we drive the quadrupole modes sufficiently strongly. We also discussed 
some of the expected experimental signatures. 

Our theory can be applied to a variety of new problems. For example, exploring the 
properties of the vortex line at finite temperatures is possible by adapting the theory we 
have used in this paper [s^. Also, the interplay between the superfiuid Mott-insulator 
transition, the quantum melting of a vortex lattice, and the quantum-Hall regime in an 
optical lattice is a very interesting topic for further research. 
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FIG. 12: This figure shows the vortex-position distribution in units of \/fi /muj.,. (a) before and (b) 
after the expansion where condensate size increased by a factor of 3.5. We used parameters J = 0.1, 
U = 100, $7 = 0, and Ng = 101. Initially vortex positions were sampled from a two-dimensional 
Gaussian distribution with an asymmetry of e = 0.9 and the product of the widths was axCy = Bq. 
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